##################################################################
# /* Author: Gautam Nair and Nicholas Sambanis */
# /* Violence Exposure and Ethnic Identification: Evidence from Kashmir */
# /* Table 19:  Power Analysis */

##################################################################
# loading packages
##################################################################
library(Hmisc)
library(foreign)
library(sandwich)
library(lmtest)
library(numDeriv)
library(stargazer)
library(ggplot2)
library(plyr)
library(gridExtra)
library(dplyr)
library(scales)

rm(list=ls())
data.working <- read.csv("nei_kashmir_replication_dataset.csv")

data.working <- as.data.frame(data.working[data.working$group1_control==1,])

dep.var <- c("identity_rank_indian_rev", 
"identity_choice_indian_both", 
"bonus_donation_india", 
"kashmir_status_india", 
"ind_pak_prefer_ind", 
"ic_index_id_no_ind_pak", 
"protest_peaceful_endorse" ,
"protest_violent_endorse")

dep.labels <-  c("Ranked Indian Identity (Reversed 1-4)", 
"Indian or Indian+Kashmiri/Kashmiri (1/0)", 
"Donation to All-India NGO (Percent)", 
"Kashmir remain Part of India (1/0)", 
"India/Pakistan: Prefer India (1/0)", 
"National Identification Index", 
"Peaceful Protest Participation (1-5)", 
"Violent Protest Endorsement (1-5)")

ind.var <- c("group_control")

# data.temp <- na.omit(data.working)

matrix.titles <- c("Control Mean", "Control SD", "SE (N=1000)", "MDE (N=1000)")

coef.matrix <- matrix(data=NA,nrow=length(dep.var),ncol=length(matrix.titles))

for(i in 1:length(dep.var)){
	
	x=0
	temp.var <- na.omit(data.working[,c(dep.var[i])])
	
	x=x+1
	mean.temp <- mean(data.working[,c(dep.var[i])], na.rm=TRUE)
	coef.matrix[i,x] <- mean.temp
	
	x=x+1
	var.temp <- var(data.working[,c(dep.var[i])], na.rm=TRUE)
	coef.matrix[i,x] <- var.temp
	
	x=x+1
	se.temp.1000 <- sqrt((var.temp/500)*2)
	coef.matrix[i,x] <- se.temp.1000
	
	x=x+1
	mde.1000 <-(0.84+1.96)*2*se.temp.1000
	coef.matrix[i,x] <- mde.1000
	
	
	my.formula <-paste(dep.var[i],'~', "religion_sunni")
	sunni.nonsunni.diff <- summary(lm(my.formula, data=data.working))$coefficient[1,2]
	x=x+1
	#coef.matrix[i,x] <- sunni.nonsunni.diff
	
}

coef.matrix <- round(coef.matrix, 2)

table <- cbind(c("Ranked Indian Identity (Reversed 1-4)", 
"Indian or Indian+Kashmiri/Kashmiri (1/0)", 
"Donation to All-India NGO (Percent)", 
"Kashmir remain Part of India (1/0)", 
"India/Pakistan: Prefer India (1/0)", 
"National Identification Index", 
"Peaceful Protest Participation (1-5)", 
"Violent Protest Endorsement (1-5)"), coef.matrix)

table <- rbind(c("Variable", "Control Mean", "Control Var.", "SE", "MDE"), table)

spectitle <- c("Power Analysis")
output.file <- c("tf_t_19_power_analysis")
columnlabels <- c("Variable", "Control Mean", "Control Variance", "SE (N=1000)", "MDE (N=1000)")
temptype <- c("latex", "text")
tempext <- c(".tex", ".txt")

for(q in 1:length(temptype)){
	temp.output <- paste(output.file, tempext[q], sep="") 
stargazer(table,
          out=temp.output,
         title= spectitle,
          font.size="footnotesize",
		summary=FALSE,
		type= temptype[q]
) 
}

